

* Open log
capture log close
log using "Data analysis\cbld-analysis01-app.log", replace text


* ************************************************************************************
* Appendix: Additional results for 'An Empirical evaluation of aggregate model output'
* ************************************************************************************

* Programme:	cbld-analysis01-app.do
* Project:		Council coalition building
* Author:		Frank Haege, Department of Politics and Administration, University of Limerick
* Contact:		frank.haege@ul.ie

* Description
*************
* This do-file generates a table and graphs presenting the results of sensitivity analyses 
* of the model's fit to the observed consensus rates in the Council. The senstitivity analyses    
* are based on models that incorporate real-world voting weights and thresholds, and in which the 
* sequence of moves of member states is ordered according to voting weights (from smallest to largest). 


* Set up Stata
version 11
clear all
macro drop _all
set linesize 80
set more off

* Generate empty matrix to save model fit statistics
matrix results = J(18,4,.)
matrix colnames results = "12" "15" "25" "Total"
matrix rownames results = Observed Null "Model 1" "PRE 1" "Model 2" "PRE 2" "Model 3" "PRE 3" /*
	*/ "Model 4" "PRE 4" "Model 5" "PRE 5" "Model 6" "PRE 6" "Model 7" "PRE 7" "Model 8" "PRE 8"
matrix list results


* Fit statistics of models with unweighted votes
************************************************

* Set starting value for matrix row number
local row = 2

* Loop through model results based on random (a) and ordered sequence of moves (b)
foreach i in "a" "b" { 

	* Load oberved consensus rates data
	use "Experiments\Other data\cbld-management05a.dta", clear
	
	* Rescale consensus rate variable
	generate oconsensus = consensus/100
	
	* Drop redundant variables
	drop consensus

	* Merge observed consensus rates data with predicted consensus rates data
	merge m:m rulemembersize using "Experiments\Experiment02\cbld-management03-exp02`i'.dta"
	drop _merge
	
	* Drop simulation results for 27 member states (under Nice and Lisbon)
	drop if year == .

	* Calculate squared deviations
	generate mdev = (oconsensus - consensus)^2
	generate ndev = (oconsensus - econsensus)^2

	* Generate sum of square
	egen mss = sum(mdev), by(nostates)
	egen tss = sum(ndev), by(nostates)
	egen tmss = sum(mdev)
	egen ttss = sum(ndev)
	drop oconsensus mdev ndev
	
	* Generate R-squared
	generate r2 = (tss - mss) / tss
	generate tr2 = (ttss - tmss) / ttss

	* Collapse data set by number of member states
		collapse (mean) r2=r2 tr2=tr2 consensus=consensus se=se lowerci=lowerci upperci=upperci /*
		*/ econsensus=econsensus ese=ese elowerci=elowerci eupperci=eupperci mss=mss tss=tss /*
		*/ tmss=tmss ttss=ttss, by(rulemembersize)

	* Label variables
	label var rulemembersize "Rule-membership period"
	label var r2 "Prop. reduction of error for regime periods"
	label var tr2 "Prop. reduction of error for entire period"
	label var consensus "Expected consensus (model)"
	label var se "Standard error (model)"
	label var lowerci "CI lower bound (model)"
	label var upperci "CI upper bound (model)"
	label var econsensus "Expected consensus (null)"
	label var ese "Standard error (null)"
	label var elowerci "CI lower bound (null)"
	label var eupperci "CI upper bound (null)"	

	* Merge with observed data
	merge 1:m rulemembersize using "Experiments\Other data\cbld-management05b.dta"
	drop _merge

	* Generate variable for running number
	generate no = _n
	replace no = no - 3.2 if no > 3

	* Merge with observed data
	merge m:m rulemembersize using "Experiments\Other data\cbld-management05c.dta"
	drop _merge

	* Update running number variable
	replace no = _n if no == .
	replace no = no - 6.1 if no > 6

	* Label running number variable
	label var no "Number of Member States"
	label def nol 1 "12" 2 "15" 3 "25", replace
	label val no nol

	* Rescale variables
	replace r2 = r2 *100
	replace tr2 = tr2 *100
	generate pconsensus = consensus *100
	generate plowerci = lowerci * 100
	generate pupperci = upperci * 100
	generate peconsensus = econsensus * 100
	generate pelowerci = elowerci * 100
	generate peupperci = eupperci * 100

	* Save R-squared measures in locals
	sort rulemembersize
	mkmat rulemembersize r2 tr2 mconsensus pconsensus peconsensus, matrix(r2)
	matrix list r2
	local r2_12 = string(round(r2[3,2], 1))
	local r2_15 = string(round(r2[6,2], 1))
	local r2_25 = string(round(r2[9,2], 1)) 
	local tr2 = string(round(r2[3,3], 1))

	* Extract results and save them in matrix
	matrix results[1,1] = r2[2,4]
	matrix results[1,2] = r2[5,4]
	matrix results[1,3] = r2[8,4]
	matrix results[2,1] = r2[3,6]
	matrix results[2,2] = r2[6,6]
	matrix results[2,3] = r2[9,6]
	local ++row
	matrix results[`row',1] = r2[3,5]
	matrix results[`row',2] = r2[6,5]
	matrix results[`row',3] = r2[9,5]
	local ++row
	matrix results[`row',1] = r2[3,2]
	matrix results[`row',2] = r2[6,2]
	matrix results[`row',3] = r2[9,2]
	matrix results[`row',4] = r2[3,3]
	matrix list results

	* Plot model and null-model estimates against observed data
	twoway (rcap pupperci plowerci no, /*
		*/ ylabel(0 20:100) xsize(6) /*
		*/ ytitle("Consensus Rate (Per Cent)") lwidth(medium) lcolor(grey) msize(medium) /*
		*/ xscale(range(0 3.5)) xlabel(1 2 3, val)) /*
		*/ (rcap peupperci pelowerci no, /*
		*/ lwidth(medium) lcolor(grey) msize(medium)) /*
		*/ (scatter peconsensus no, msymbol(s) mcolor(white)) /*
		*/ (scatter pconsensus no, msymbol(d) mcolor(white)) /*
		*/ (scatter mconsensus no, msymbol(dh) mcolor(white)) /*
		*/ (scatter consensus1994 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus1995 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus1996 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus1997 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus1998 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus1999 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2000 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2001 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2002 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2003 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2004 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2005 no, msymbol(oh) mcolor(white)) /*
		*/ (scatter consensus2006 no, msymbol(oh) mcolor(white)) /*
		*/ , legend(order(- "Predicted:" 4 "Computational Model" 3 "Null Model" /*
		*/ - " " - "Observed:" 5 "Mean" 6 "Individual Values" - " " - "Total {it:PRE}: `tr2'%") /*
		*/ cols(1) ring(1) position(3) region(lcolor(black))) scale(1.2) /*
		*/ text(10 .3 "{it:PRE} :") text(10 1 "`r2_12'%") text(10 2 "`r2_15'%") text(10 3 "`r2_25'%") /*
		*/ saving("Data analysis\Graphs\cbld-analysis01-app01`i'.gph", replace)
}


* Fit statistics of models with weighted votes and several voting thresholds
****************************************************************************

* Loop through model results based on 1 to 3 voting thresholds
foreach j of num 1 2 3 { 
	* Loop through model results based on random (a) and ordered sequence of moves (b)
	foreach i in "a" "b" {

		* Load observed consensus rates data
		use "Experiments\Other data\cbld-management05a.dta", clear
		
		* Rescale consensus rate variable
		generate oconsensus = consensus/100
		
		* Drop redundant variables
		drop consensus

		* Merge observed consensus rates data with predicted consensus rates data
		merge m:m rulemembersize using "Experiments\Experiment03\cbld-management04-exp03`i'`j'.dta"
		drop _merge
	
		* Drop simulation results for 27 member states (under Nice and Lisbon)
		drop if year == .
	
		* Calculate squared deviations
		generate mdev = (oconsensus - consensus)^2
		generate ndev = (oconsensus - econsensus)^2

		* Generate sum of square
		egen mss = sum(mdev), by(nostates)
		egen tss = sum(ndev), by(nostates)
		egen tmss = sum(mdev)
		egen ttss = sum(ndev)
		drop oconsensus mdev ndev

		* Generate R-squared
		generate r2 = (tss - mss) / tss
		generate tr2 = (ttss - tmss) / ttss

		* Collapse dataset by number of member stats
		collapse (mean) r2=r2 tr2=tr2 consensus=consensus se=se lowerci=lowerci upperci=upperci /*
			*/ econsensus=econsensus ese=ese elowerci=elowerci eupperci=eupperci mss=mss tss=tss /*
			*/ tmss=tmss ttss=ttss, by(rulemembersize)

		* Label variables
		label var rulemembersize "Rule-membership period"
		label var r2 "Prop. reduction of error for regime periods"
		label var tr2 "Prop. reduction of error for entire period"
		label var consensus "Expected consensus (model)"
		label var se "Standard error (model)"
		label var lowerci "CI lower bound (model)"
		label var upperci "CI upper bound (model)"
		label var econsensus "Expected consensus (null)"
		label var ese "Standard error (null)"
		label var elowerci "CI lower bound (null)"
		label var eupperci "CI upper bound (null)"	

		* Merge with observed consensus rates data
		merge 1:m rulemembersize using "Experiments\Other data\cbld-management05b.dta"
		drop _merge

		* Generate variable for running number
		generate no = _n
		replace no = no - 3.2 if no > 3

		* Merge with observed consensus rates data
		merge m:m rulemembersize using "Experiments\Other data\cbld-management05c.dta"
		drop _merge

		* Update running number variable
		replace no = _n if no == .
		replace no = no - 6.1 if no >= 7

		* Label running number variable
		label var no "Number of Member States"
		label def nol 1 "12" 2 "15" 3 "25", replace
		label val no nol

		* Rescale variables
		replace r2 = r2 *100
		replace tr2 = tr2 *100
		generate pconsensus = consensus *100
		generate plowerci = lowerci * 100
		generate pupperci = upperci * 100
		generate peconsensus = econsensus * 100
		generate pelowerci = elowerci * 100
		generate peupperci = eupperci * 100

		* Save R-squared statistics in locals
		sort rulemembersize
		mkmat rulemembersize r2 tr2 mconsensus pconsensus peconsensus, matrix(r2)
		matrix list r2
		local r2_12 = string(round(r2[3,2], 1))
		local r2_15 = string(round(r2[6,2], 1))
		local r2_25 = string(round(r2[9,2], 1)) 
		local tr2 = string(round(r2[3,3], 1))

		* Extract results and save them in matrix
		matrix results[1,1] = r2[2,4]
		matrix results[1,2] = r2[5,4]
		matrix results[1,3] = r2[8,4]
		matrix results[2,1] = r2[3,6]
		matrix results[2,2] = r2[6,6]
		matrix results[2,3] = r2[9,6]
		local ++row
		matrix results[`row',1] = r2[3,5]
		matrix results[`row',2] = r2[6,5]
		matrix results[`row',3] = r2[9,5]
		local ++row
		matrix results[`row',1] = r2[3,2]
		matrix results[`row',2] = r2[6,2]
		matrix results[`row',3] = r2[9,2]
		matrix results[`row',4] = r2[3,3]
		matrix list results

		* Plot model and null-model estimates against observed data
		twoway (rcap pupperci plowerci no, /*
			*/ ylabel(0 20:100) xsize(6) /*
			*/ ytitle("Consensus Rate (Per Cent)") lwidth(medium) lcolor(grey) msize(medium) /*
			*/ xscale(range(0 3.5)) xlabel(1 2 3, val)) /*
			*/ (rcap peupperci pelowerci no, /*
			*/ lwidth(medium) lcolor(grey) msize(medium)) /*
			*/ (scatter peconsensus no, msymbol(s) mcolor(white)) /*
			*/ (scatter pconsensus no, msymbol(d) mcolor(white)) /*
			*/ (scatter mconsensus no, msymbol(dh) mcolor(white)) /*
			*/ (scatter consensus1994 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus1995 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus1996 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus1997 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus1998 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus1999 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2000 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2001 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2002 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2003 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2004 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2005 no, msymbol(oh) mcolor(white)) /*
			*/ (scatter consensus2006 no, msymbol(oh) mcolor(white)) /*
			*/ , legend(order(- "Predicted:" 4 "Computational Model" 3 "Null Model" /*
			*/ - " " - "Observed:" 5 "Mean" 6 "Individual Values" - " " - "Total {it:PRE}: `tr2'%") /*
			*/ cols(1) ring(1) position(3) region(lcolor(black))) scale(1.2) /*
			*/ text(10 .3 "{it:PRE} :") text(10 1 "`r2_12'%") text(10 2 "`r2_15'%") text(10 3 "`r2_25'%") /*
			*/ saving("Data analysis\Graphs\cbld-analysis01-app02`i'`j'.gph", replace)
	}
}

* TABLE A1.	Observed vs. predicted consensus rates for various model specifications
***********************************************************************************
matrix list results, format(%2.0f)


log close
exit


* Fig. A1. Unweighted votes, random sequence of moves (base model)
******************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app01a.gph"

* Fig. A2. Unweighted votes, sequence of moves ordered by size
**************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app01b.gph"

* Fig. A3. Weighted votes, random sequence of moves, one threshold under Nice
*****************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02a1.gph"

* Fig. A4. Weighted votes, sequence of moves ordered by size, one threshold under Nice
**************************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02b1.gph"

* Fig. A5. Weighted votes, random sequence of moves, two thresholds under Nice
******************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02a2.gph"

* Fig. A6. Weighted votes, sequence of moves ordered by size, two thresholds under Nice
***************************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02b2.gph"

* Fig. A7. Weighted votes, random sequence of moves, three thresholds under Nice
********************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02a3.gph"

* Fig. A8. Weighted votes, sequence of moves ordered by size, three thresholds under Nice
*****************************************************************************************
graph use "Data analysis\Graphs\cbld-analysis01-app02b3.gph"
